Microbial species pool-mediated diazotrophic community assembly in crop microbiomes during plant development

ABSTRACT Plant-associated diazotrophs strongly relate to plant nitrogen (N) supply and growth. However, our knowledge of diazotrophic community assembly and microbial N metabolism in plant microbiomes is largely limited. Here we examined the assembly and temporal dynamics of diazotrophic communities across multiple compartments (soils, epiphytic and endophytic niches of root and leaf, and grain) of three cereal crops (maize, wheat, and barley) and identified the potential N-cycling pathways in phylloplane microbiomes. Our results demonstrated that the microbial species pool, influenced by site-specific environmental factors (e.g., edaphic factors), had a stronger effect than host selection (i.e., plant species and developmental stage) in shaping diazotrophic communities across the soil–plant continuum. Crop diazotrophic communities were dominated by a few taxa (~0.7% of diazotrophic phylotypes) which were mainly affiliated with Methylobacterium, Azospirillum, Bradyrhizobium, and Rhizobium. Furthermore, eight dominant taxa belonging to Azospirillum and Methylobacterium were identified as keystone diazotrophic taxa for three crops and were potentially associated with microbial network stability and crop yields. Metagenomic binning recovered 58 metagenome-assembled genomes (MAGs) from the phylloplane, and the majority of them were identified as novel species (37 MAGs) and harbored genes potentially related to multiple N metabolism processes (e.g., nitrate reduction). Notably, for the first time, a high-quality MAG harboring genes involved in the complete denitrification process was recovered in the phylloplane and showed high identity to Pseudomonas mendocina. Overall, these findings significantly expand our understanding of ecological drivers of crop diazotrophs and provide new insights into the potential microbial N metabolism in the phyllosphere. IMPORTANCE Plants harbor diverse nitrogen-fixing microorganisms (i.e., diazotrophic communities) in both belowground and aboveground tissues, which play a vital role in plant nitrogen supply and growth promotion. Understanding the assembly and temporal dynamics of crop diazotrophic communities is a prerequisite for harnessing them to promote plant growth. In this study, we show that the site-specific microbial species pool largely shapes the structure of diazotrophic communities in the leaves and roots of three cereal crops. We further identify keystone diazotrophic taxa in crop microbiomes and characterize potential microbial N metabolism pathways in the phyllosphere, which provides essential information for developing microbiome-based tools in future sustainable agricultural production.

P lants live in close association with diverse microorganisms (i.e., the plant microbiota), which are indispensable for plant fitness and performance in many aspects (1)(2)(3)(4).Among them, symbiotic and non-symbiotic diazotrophs (nitrogen fixers, e.g., Rhizobium and Azospirillum) not only significantly contribute to plant nitrogen (N) supply but also benefit plant growth promotion via provisions of primary and secondary metabolites (1,5,6).Hence, diazotrophs are considered ideal candidates for sustaining sustainable agricultural production as plant growth-promoting rhizobacteria (PGPR) (6,7).Increasing evidence showed that biological N-fixation (BNF) in roots, stalks, and even plant leaf surfaces and endospheres contributes significantly to the N nutrition for the host plant and the N budget in terrestrial ecosystems (5,(8)(9)(10)(11)(12).However, current knowledge about the diversity and ecological drivers of diazotrophs is mainly derived from soils.Given the processes of the plant microbiome assembly differ from those of the soil microbiome (e.g., the effect of host factors), it is critical to identify the fundamental ecological processes that govern the distribution and dynamics of diazotrophic communities across the soil to plant habitats.This knowledge could provide important information for harnessing the plant microbiome for sustainable agricultural production (7,13).
Recent studies have demonstrated the immense role of multiple host and environ mental factors in driving the assembly of plant-associated microbiomes (14)(15)(16)(17)(18).For example, it was suggested that the assembly of both bacterial and fungal communities across the soil-plant continuum is predominantly determined by plant compartment and host species, followed by environmental factors (19,20).As for a specific compart ment like the phylloplane (leaf surface) and rhizoplane, microbiomes are co-influenced by the host (e.g., plant species) and environmental (e.g., site, soil properties, and climate) factors (15,(21)(22)(23)(24).These studies, however, mainly focused on overall bacteria and fungi, and a significant knowledge gap exists on how these multiple host and environmental factors shape the diversity and structure of diazotrophic communities along the soilplant continuum.Given the importance of diazotrophs in plant N supply and growth promotion, theoretically, a close relationship between host and associated diazotrophs is expected (e.g., plant-specific keystone taxa) but empirical evidence is lacking.Numer ous studies on soil and rhizosphere samples have suggested that site, soil pH, and fertilization are the main driving factors determining the distribution of diazotrophic communities in agriculture ecosystems (25)(26)(27)(28)(29)(30)(31).It was suggested that long-term N fertilization drastically reduced N fixation and suppressed certain diazotrophic taxa in soils (32).However, the influence of agricultural management regimes on crop-associ ated diazotrophs remains largely unknown.
In addition to their role in fixing N 2 , emerging evidence is showing that microbiomes might play a role in nitrification, denitrification, and N assimilation, and contribute to the decomposition of methanol and isoprene in the phyllosphere (i.e., both epiphytic and endophytic microhabitats of plant aerial parts) (15,(33)(34)(35).For example, metagenomic analyses at the contig level showed that the maize phyllosphere microbiome harbors diverse N-cycling functional genes (15).Phyllosphere microorganisms are expected to play critical roles in atmospheric N cycling since the phyllosphere acts as the plantatmosphere interface and is usually exposed to multiple environmental changes like N deposition (36,37).However, a systematic understanding and genome-based evidence on the role of phyllosphere microbiomes in N cycling in agroecosystems is still scarce.
In this study, we used amplicon sequencing of the nifH gene to examine diazo trophic communities across soils (bulk and rhizosphere soils), plastic leaf (represent ing local air environment), and multiple compartments (rhizoplane, root endosphere, phylloplane, and leaf endosphere) of three crops (maize, wheat, and barley) under four fertilization regimes in two major agricultural areas of China (Fig. S1).We then explored the N-cycling pathways of phylloplane microbiomes via metagenomic binning.We aimed to (1) determine the distribution pattern of diazotrophic communities along the soil-plant continuum of different crops and across different plant developmental stages (2); identify dominant diazotrophic taxa of crop microbiomes and their poten tial ecological role in microbial co-occurrence network and plant performance; and (3) examine the potential microbial N metabolism pathways in phylloplane microbiomes in agricultural ecosystem, based on metagenome-assembled genomes.We hypothesized that (1) both host selection, driven by crop species and developmental stage, and site-specific microbial species pool will influence equally the assembly of crop-associated diazotrophic communities (2); some specific diazotrophic taxa will be associated with plant selection effects and may play key roles in microbial co-occurrence network and crop production; and (3) phylloplane microbiomes may possess ample N metabolism pathways, as they are typically exposed to atmospheric N cycling.
We then explored the temporal dynamics of the diazotrophic community.The diazotrophic community in the phylloplane (R 2 = 16.8%,P < 0.001) was more sensitive to the plant developmental stage than those in the rhizoplane (R 2 = 6.8%,P = 0.03) and rhizosphere (R 2 = 2.5%, P = 0.015) (Fig. 1c; Table S1).Sampling time (R 2 = 36.4%,P < 0.001) also showed a significant effect on the diazotrophic community in the plastic leaf (Fig. 1c; Table S1).For the alpha-diversity, the diazotrophic Shannon index was highest in bulk and rhizosphere soils, and followed by rhizoplane, but showed no significant difference among the three crops (Fig. S2).Sampling time had a significant effect on the Shannon index of diazotrophic communities in both air and rhizoplane samples (Fig. S2).
To further explore the response of diazotrophic communities in each compart ment to four fertilization practices, that is, no N fertilizer (Control), 80% N fertilizer (80%N), biochar treatment (80%NSB), and straw treatment (80%NS), sub-dataset of each compartment were subjected to PERMANOVA analyses.It showed that fertilization practice explained 33%-51% of community variation in each compartment at two sites (Fig. S3a and b).Particularly, biochar treatment displayed a significant effect on diazotrophic community composition in most plant and soil compartments (Fig. S3a and  b), and significantly increased diazotrophic Shannon index in both rhizosphere and bulk soils at the site XC, compared with its counterpart treatment (80%NS) (Fig. S2e).

Taxonomic composition and dominant taxa of crop-associated diazotrophs
Taxonomic classification showed that diazotrophic community largely varied between plant and soil habitats, with more abundant Methylobacterium in plant compartments (leaves, roots, and grain) (~15.3%vs ~2.0% in soils) while more abundant Bradyrhizobium in soils (~44.3%vs ~18.7% in plants) (Fig. 2b; Fig. S3c).As for two sites (including soil and air samples), site XC had a higher proportion of Azospirillum and Sinorhizobium while site QJ harbored more abundant Methylobacterium, Methylosoma, and Rhodopseudomonas (Fig. 2c; Fig. S3d).Among the three crops (only includes plant samples), more abundant Rhizobium was observed in barley while more abundant Pseudacidovorax was found in maize or wheat (Fig. 2d; Fig. S3e).The linear mixed model analysis suggested that the plant developmental stage had the strongest effect on Rhizobium (F = 9.7, P = 9.7e-5), dissimilarity was tested with PERMANOVA."Niche" represents the effect of compartment niche, "Crop" represents the effect of crop species, "Site" represents the effect of site, and "Stage" represents the effect of developmental stage."XC" represents the site "Xuchang, " "QJ" represents the site "Qujing." which was more abundant at the tasseling stage than at other stages in both plastic leaf and maize phylloplane (Fig. 2e; Fig. S3f ).The differential abundance analysis showed that some members within Bradyrhizobium were significantly enriched in 80%NSB treatment in soils in both XC (OTU9001 and OTU19639) and QJ (OTU6065) sites (Fig. S4).
To further identify the keystone diazotrophic taxa across the soil-plant continuum, we defined the dominant taxa as the OTUs with a relative abundance ≥0.1% and present in at least one-third of samples evaluated and explored the composition of dominant sub-communities.It showed that crop diazotrophic communities were dominated by a few taxa (0.4%-0.9% of phylotypes) mainly belonging to Methylobacterium, Azospirillum, Bradyrhizobium, and Rhizobium, which accounted for 57%-72% of the total sequences in each crop (Fig. 3a and b).By contrast, soil diazotrophic communities were dominated by a few taxa mainly belonging to Bradyrhizobium and Azospirillum (Fig. 3b).Further analysis showed that eight dominant taxa belonging to Azospirillum (7 OTUs) and Methylobac terium (1 OTU) were commonly present in the three crops with close phylogenetic relationships and high relative abundance (17%-27%) in each crop, but were absent in soils, and thus were defined as crop-specific taxa (termed the A&M clade) (Fig. 3a).Mantel test based on maize season samples indicated that soil pH, NO 3 --N, and DOC were important drivers of the diazotrophic community and that maize dominant sub-com munity in the rhizoplane and rhizosphere had greater correlations with rhizosphere nitrogenase activity than whole and non-dominant sub-communities across three plant developmental stages (Table S2).Random forest modeling analyses indicated that the relative abundance of the A&M clade taxa was the most important predictor for crop yield in each crop species (Fig. S5a through c).

Co-occurrence networks of crop-associated diazotrophic communities
Network analysis showed that the co-occurrence patterns of the diazotrophic commun ity differed among the three crops, with more nodes belonging to Azospirillum and Sinorhizobium in maize and wheat, and Bradyrhizobium and Rhizobium in barley (Fig. 4).More importantly, all eight A&M clade taxa had significantly higher network degree than any other dominant and non-dominant taxa (Fig. 4a and c).Among the top 10 hubs with the highest degree in the network of each crop, 7-8 hubs belonged to the A&M clade (Fig. 4d), among which OTU17082 (Methylobacterium) was identified as the most important network hubs for both maize and wheat (Fig. 4b and d).Interestingly, the hubs belonging to A&M clade taxa showed a positive correlation between them but mainly showed a negative correlation with other diazotrophic taxa in the network of three crops (Fig. 4b and d).In addition, diazotrophic co-occurrence patterns differed among different compartments, with more nodes belonging to Azospirillum in the phylloplane (both plant and plastic) (Fig. S6).

Compartment and site dominate over host species in determining the diazotrophic community across the soil-plant continuum
Understanding the ecological drivers of crop-associated diazotrophic communities is integral to harness diazotrophs for sustainable agricultural production.To reach this target, we examined diazotrophic communities across the soil-plant continuum of three crops under different fertilization practices in two geographically distant soils and incorporated artificial plants as background controls.Our study demonstrated that compartment and site are the primary drivers shaping the assembly of crop diazotrophic communities across the soil-plant continuum.As for the bulk and rhizosphere soils, diazotrophic communities were largely determined by site (77%-85%).These findings were consistent with what has been observed previously in bacterial and fungal communities (19,20,38,39).Surprisingly, we found that the site (10%-62%) plays a more important role in determining diazotrophic community in most plant compartments than crop species (8%-21%) and plant developmental stage (6%-16%).More specifically, maize diazotrophs in the root endosphere from site QJ were more closely assembled with barley diazotrophs from the same site, rather than with same the crop from site XC (Fig. 1b).This result contrasts with our previous observation showing that crop species (bacteria ~52%, fungi ~41%) dominated over site (bacteria ~20%, fungi ~15%) in determining the assembly of bacterial and fungal communities in plant compartments (19,20).The diazotrophic community in the artificial leaf (i.e., aerobic background) differed between the two sites and was more sensitive to sampling time than those of the maize, indicating the important role of season-dependent climatic factors on leaf diazotrophs.These results collectively reveal that, in comparison to the overall bacterial and fungal communities, the diazotrophic community, as a functional group within bacteria, is more sensitive to site-dependent environmental factors but more stable to host effect caused by crop species and developmental stage.
Thus, we proposed a conceptual model to explain the ecological mechanisms that govern the assembly of crop-associated diazotrophic communities under different environmental conditions (Fig. 6).Our study demonstrated that the influence of site on diazotrophic communities outweighs that of crop species and developmental stage in bulk/rhizosphere soils and most plant compartments.Further analyses indicated that soil and airborne diazotrophs were important sources of plant diazotrophs.Together, these results revealed that the local microbial species pool developed under strong site effects plays an important role in shaping diazotrophic communities across the soil-plant continuum.Our results showed that soil pH and nutrients including availa ble N and DOC are primary determinants affecting the assembly of soil diazotrophic communities.Site-dependent climatic factors (e.g., temperature and precipitation) and agronomic regimes (e.g., fertilization habits) also play key roles in structuring the microbial species pool of diazotrophs in soils and atmosphere, which further exerts influence on the diazotrophic communities in plant compartments (Fig. 6).On the other hand, the site-dependent edaphic and climatic factors also define plant growth and development (2,40), which might, in turn, shape crop diazotrophic community assembly.The sites QJ and XC had distinct climates (subtropical monsoon vs temperate mon soon), soil properties (acidic red soil vs alkaline Fluvo-aquic soil), and different atmos pheric deposition (25,41), which facilitate the differentiation of microbial species pool between two sites and profoundly affected the assembly of crop-associated diazotrophic communities (Fig. 6).
Furthermore, many studies have reported that different agricultural management regimes (e.g., fertilization practice) can significantly affect the diversity and structure of soil diazotrophic communities (29,32,42).Here we provide new evidence that soil amendment practices like biochar addition can not only influence the diazotro phic community in soils (enhancing diversity, changing the community structures, and enriching Bradyrhizobium) but also affect the diazotrophic community in plant compartments.These results concur with previous findings on the positive impacts of biochar on soil microbiomes and fertility (43,44) and extend our understanding of crop microbiomes.The biochar amendment may exert its influence on crop diazotrophs directly by affecting the soil microbial species pool, or indirectly by affecting plant growth and development (45).These findings implied the possibility of manipulating the diazotrophic community for sustainable crop production and thus has significant implications for agricultural practices.

Keystone diazotrophic taxa in soil and crop microbiomes
Our study further showed that the crop diazotrophic community was dominated by a few phylotypes which were mainly affiliated with Methylobacterium, Azospirillum, Bradyrhizobium, and Rhizobium.Coincidently, previous studies suggest that bacterial and fungal communities in soils (46) and plant compartments (19,20) were dominated by a few dominant taxa.More intriguingly, eight dominant taxa (A&M clade) comprised of Azospirillum and Methylobacterium were identified as keystone diazotrophic taxa in three crops and played a potential role in predicting crop yield in this study.Similarly, a recent study suggested that the diversity of keystone phylotypes in soil microbial communities had positive associations with crop production (27).We further found that the network hubs belonging to the A&M clade had higher proportions of positive edges between them (i.e., positive associations).Coincidently, a previous study by Madhaiyan et al. found that the co-inoculation of Methylobacterium oryzae and Azospirillum brasilense displayed a significantly positive effect on the growth and nutrient uptake of various crops, including tomato, red pepper, and rice (47).These together with our findings suggested that Azospirillum and Methylobacterium might synergistically function in plant compartments.However, we must acknowledge that the positive correlations between A&M clade taxa observed in the network may be attributable to similar selection factors, rather than actual synergistic interactions.The real microbial interactions within the A&M clade and among different groups need further investigation.The observed diversity and dominance of diazotrophic communities across various crops further revealed the functional redundancy of plant microbiomes, which might be a strategy of the micro biome to warrant a given function being maintained in a changing environment.
Furthermore, these A&M clade taxa were significantly enriched in plant compart ments in comparison to bulk and rhizosphere soils, suggesting their preference for plant compartments.Some members of Methylobacterium were frequently identified as dominant taxa in the phyllosphere of various plants (15,48) and demonstrated to be able to fix N 2 and promote plant growth (49).Our metagenomic analyses further showed that two MAGs belonging to Methylobacterium in the phylloplane contained genes involved in N and methane metabolisms, further confirming their potential function in the phyllosphere.The members of Azospirillum are capable of colonizing a wide range of habitats like the phyllosphere and rhizosphere and play a key role in modulating host fitness (50,51), which may explain their dominant distribution in various plant compart ments.These results suggest that Methylobacterium and Azospirillum may play vital ecological roles in plant performance while their true function needs further validation through synthetic community experiments.
We further found that Rhizobium was more abundant in barley than in maize or wheat, suggesting that the barley species may preferentially enrich members of Rhizobium.Except for N fixation with leguminous plants, many members of Rhizobium are supposed to be important PGPR which cannot form nodules with plants but are capable of promoting plant growth via secreting phytohormones and siderophores, inducing plant immunity response and fitness to biotic and abiotic stresses on nonlegume plants (52)(53)(54).As for soils, dominant diazotrophic sub-communities were dominated by Bradyrhizobium, which was further enriched by biochar amendment (80%NSB).Bradyrhizobium was consistently recognized as the dominant diazotrophic taxa across diverse ecosystems and is believed to play key roles in active N fixation in soils (27,55,56).The identification of these keystone diazotrophic taxa provides critical information for engineering synthetic communities (SynComs) and manipulating cropdiazotrophs interactions for future sustainable agricultural production.On the other hand, we must acknowledge that the identification of these dominant and keystone taxa might be influenced by the bias of primers used in this study (57)(58)(59), as previous studies have indicated that the use of primer PolF/PolR for nifH gene amplification may result in a lower diversity and a preference for amplifying Proteobacteria, Firmicutes, and Actinobacteria (57,59).More intensive research employing different primers is necessary to better profile soil-and plant-associated diazotrophic communities in the future.

The potential microbial N metabolisms and function in the phylloplane
Phyllosphere microbiomes can provide benefits to the plant in multiple aspects (60)(61)(62), but our knowledge of microbial N metabolism in the phyllosphere is largely limited.Our metagenomic binning for phylloplane microbiomes suggested that 48 of 58 retrieved MAGs harbored genes related to electron transfer (e.g., nifF/J) and cofactor (e.g., nifU/V/B) in N fixation processes.However, none of the 58 MAGs possessed genes encoding nitrogenase (i.e., nifH/D/K), and some dominant diazotrophic taxa identified by the amplicon sequencing, such as Azospirillum, Bradyrhizobium, and Rhizobium, were not recovered in metagenomic binning (i.e., not within 58 MAGs).These might have been caused by the limited sequencing depth in this study and the low MAG recov ery rate due to the low abundance of diazotrophs in the phylloplane.The activity and N fixation function of the diazotrophic community in the phylloplane therefore need further exploration based on deeper sequencing analysis or more powerful techniques such as 15 N tracing or stable isotope probe (SIP).Remarkably, a high-qual ity MAG 25.A with high identity similarity to P. mendocina possessed genes involved in a complete denitrification pathway (i.e., covering all four denitrification processes from NO 3 − →NO 2 − →NO→N 2 O→N 2 ).Though biological nitrification and N 2 O production via ammonia-oxidizing bacteria (AOB)-mediated nitrification in the phyllosphere have been previously suggested in forests (34) and pasture (63), to the best of our knowl edge, this is the first time that microbial genomes with complete denitrification-related genes were identified in the phylloplane.Although the genome contains genes indicating potential for certain functions, this does not necessarily mean that the microbe expresses these functions in the phylloplane.This observation suggested that phylloplane microbiomes possess diverse metabolic pathways that confer their plasticity to adapt to different environments and conditions.Their actual ecological function in the phylloplane requires further experimental validation.On the other hand, the potential denitrification at the site QJ could be explained by the fact that excessive atmospheric N deposition in the form of NOx (34,64).N deposition is one of the major N sources for agricultural areas in southwest China close to the QJ site (~20 kg N ha −1 yr −1 ) (64,65).High denitrification ability under aerobic condition had been reported in P. mendocina from sewage treatment and other members of Pseudomonas (66)(67)(68).These findings, together with the detection of genes related to aerobic respiration, fumarate reductase, lactate dehydrogenase, formate dehydrogenase, and non-ribosomal peptide synthetase (NRPS) in MAG 25.A, suggest that P. mendocina has the potential for diverse metabolic pathways.
We further identified a complete MAG 13.A that was highly identical to Providencia spp.and possessed genes related to complete DNRA pathways and photosystem.The intermediate product of denitrification and DNRA processes, nitric oxide (NO), is an important gaseous signaling molecule that regulates multiple physiological processes associated with plant development and stress responses, such as stomatal closure and photosynthesis regulation (69,70).Therefore, nitrate reduction-related functional genes and bacterial taxa retrieved in this study might be closely associated with plant signaling in the phyllosphere, which deserves exploration in further study.Overall, this study significantly advances our knowledge of the microbial N metabolism and the potential ecological function of phyllosphere microbiomes.

Conclusions
This study provides the first comprehensive empirical understanding of the diversity, composition, and keystone taxa of diazotrophic communities across multiple soil and plant compartments of three cereal crops.Our results demonstrated that the assembly of diazotrophic communities across the soil-plant continuum is determined largely by compartment and site, followed by crop species and developmental stage.We further found that the effect of the site on diazotrophic communities outweighed those of crop species and developmental stage in most of the soil and plant compartments and that soil and airborne diazotrophs were important sources of plant diazotrophs.These findings highlight that the site-specific microbial species pool developed under long-term site effect plays a primary role in shaping the diazotrophic communities across the soil-plant continuum.Our study further identified eight dominant diazotrophic taxa (belonging to Azospirillum and Methylobacterium) as keystone taxa in crop diazotrophs, which might be responsible for regulating microbial network stability and acting as key N fixers and plant growth-promoting bacteria in crops.
Furthermore, 58 bacterial MAGs were recovered from plastic and maize phylloplane microbiomes and the majority of them were identified as novel species (37 MAGs) and harbored genes related to N metabolism (e.g., nitrate reduction) processes.We further identified two high-quality MAGs (close to P. mendocina and Providencia spp.) harboring genes involved in the complete denitrification pathway and DNRA pathway, respectively.These findings significantly advanced our knowledge of the assembly of diazotrophic communities across soil and plant compartments and provided an important direction for harnessing the plant microbiome for sustainable development goals.

Experimental design and sampling
The field experiments with seven fertilization treatments were established in Xuchang, Henan province (XC, the alkaline Fluvo-aquic soil located in northern China), and Qujing, Yunnan province (QJ, the acidic red soil located in southwest China) in 2016 as previously described (19,20).To evaluate the distribution pattern of diazotrophic communities along the soil-plant continuum of different crops and across different plant develop mental stages, four treatments including (T1) Control (zero N fertilizer); (T2) 80%N (20% N reduction based on local farmers' N rate); (T3) 80%NS (80%N treatment plus straw covering); and (T4) 80%NSB (80%NS treatment plus biochar addition) were involved in this study.For each crop, samples from five compartments (including leaf, grain, root, rhizosphere, and bulk soil) were collected following our previous method (15,19).For maize season, sample collections were performed in June, August, and September of 2017, corresponding to maize seedling, tasseling, and mature stage, respectively (n = 282).Sampling for the wheat/barley season was performed once at the grain filling stage in April 2018 (n = 162) (Fig. S1a and b).Artificial plants made of plastic material were planted as "background controls" in the field, and plastic leaves were sampled at the same time with soil and plant sample collection (n = 24) (15).Thus, a total of 468 samples were retrieved for this study, of which maize season samples in 2017 were used to examine bacterial and fungal communities in our previous studies (15,19,20).Soil physicochemical parameters (e.g., pH, NH 4 + -N, and NO 3 --N) were measured according to previous protocols (15,20,71,72).Soil nitrogenase activity was measured using the acetylene reduction method (15,20).The crop yield data were collected based on the field harvest (15).

DNA extraction and nifH gene amplification
Rhizosphere and bulk soil DNA were extracted from 0.4 g soil using the PowerSoil DNA Isolation Kit (MO BIO Laboratories, Carlsbad, CA, USA).Epiphytic and endophytic DNA from leaves (10-15 g) and roots (3-5 g) of maize and wheat/barley, and microbial DNA from the plastic leaf and plant grain samples were extracted using the PowerSoil DNA Isolation Kit as described previously (15,19).The nifH genes were amplified using the primer pairs PolF/PolR (41,73).In total, 426 DNA samples from soils, multiple plant compartments, and plastic leaves were successfully amplified (42 DNA samples failed to be amplified, Fig. S1).More detailed descriptions on the field trial and the nifH gene amplification are available in Supplementary Information (Method S1).

Bioinformatic analysis for amplicon sequencing
Raw sequences were analyzed using USEARCH (v10.0)(74) and QIIME v1.91 (75) as previously described (19,32).Briefly, paired-end sequences were merged into a single sequence and the primer sequences were trimmed, and the resultant sequences were quality-filtered (maximum expected error 0.5) in USEARCH.The high-quality sequences were aligned against the nifH gene database to remove chimeras by UCHIME (76).OTU clustering was performed at 95% similarity using UCLUST and all singletons were removed (74).The sequences were further converted to amino acid sequences using the FunGene Pipeline of the Ribosomal Database Project, and the amino acid sequen ces that did not match the nifH protein sequence or contained termination codons were removed (56).The remaining representative sequences were annotated using the nifH gene database (77).This paired-end sequencing approach resulted in 17,773,861 high-quality reads, which were assembled into 18,750 OTUs.The OTU table was rarefied to 1,380 reads per sample for alpha-diversity estimates, and the cumulative sum scaling (CSS) method was for normalization before beta-diversity analyses.

Statistical analysis
A nonparametric statistical test (Kruskal-Wallis test or Wilcoxon test) was used to evaluate the alpha-diversity and the taxonomical difference of the diazotrophic community observed under different host and environmental factors.The beta-diversity was assessed by computing weighted UniFrac distance matrices and visualized using non-metric multi-dimensional (NMDS) ordinations.The relative contribution of host and environmental factors on community dissimilarity was tested with PERMANOVA or nested PERMANOVA using the adonis function of the "vegan" R package (91) based on weighted UniFrac distances.Source-tracking analyses were carried out to determine the potential sources of crop diazotrophs using SourceTracker (v.1.0)(92), following the methodology described in previous studies (15,19).The phylogenetic tree of dominant taxa was built by FastTree and visualized in iTOL (87).Differential abundance analysis was performed using EdgeR's generalized linear model (GLM) approach (93).This analysis reveals which diazotrophic taxa were significantly enriched in different compartments or under various treatments.Mantel test was performed to explore Spearman's correla tions between diazotrophic communities, soil physicochemical characteristics, and soil enzyme activities using the "vegan" package (20).The network analysis was performed using CoNet (94) in Cytoscape v3.5 (95), and only robust (Spearman's r > 0.6 or r < -0.6) and statistically significant (P < 0.01) correlations were retained for visualization in Gephi (96).Random forest analyses were performed to identify the relative importance of different sub-communities in predicting crop yield using the "randomForest" R package (15,97).All statistical analyses were carried out in R (http://www.r-project.org).More information on the PERMANOVA, network analysis, and random forest modeling analysis is detailed in our previous publication (15,19).

FIG 1
FIG 1 Distribution and temporal dynamics of crop-associated diazotrophic communities.(a) Non-metric multi-dimensional (NMDS) ordinations depicting the distribution patterns of diazotrophic communities along the soil-plant continuum from all samples (n = 426), maize samples (n = 270), and wheat/barley samples (n = 156).Ellipses include 80% of samples from each site.(b) NMDS ordinations showing the effects of crop species and site on diazotrophic communities in each compartment (phylloplane: n = 60; other niches: n = 48).The gradient bar chart represents the influence of the site on the diazotrophic communities, with intensified red indicating an increased influence of the site on the community.(c) NMDS ordinations of maize diazotrophic communities across three plant developmental stages in the phylloplane, rhizoplane, and rhizosphere (phylloplane: n = 144; other niches: n = 54).The relative contribution of host and environmental factors on community

FIG 3
FIG 3 Taxonomic composition and distribution patterns of dominant sub-community in soils and plant compartments.(a) Identification of dominant taxa for three crop hosts (only considering leaf, root, and grain compartments; maize: n = 168, wheat: n = 45, barley: n = 45) and soils (include rhizosphere and bulk soils, n = 132).The numbers above the A&M clade taxa represent its average relative abundance in plant compartments.(b) Pie plot showing the proportion and community composition of the dominant sub-community in each crop and soil.The dominant taxa of diazotrophic communities were defined as OTUs with an average relative abundance ≥0.1% and present in at least one-third of the samples evaluated.The genera with the relative abundance <1% are grouped into "Other."

FIG 5
FIG 5 Phylogenetic tree and nitrogen-cycling pathways of the MAGs recovered from phylloplane microbiomes.(a) Phylogenetic tree of 58 bacterial MAGs recovered from phylloplane microbiomes, including 31 MAGs recovered from the plastic leaf (named ".A") and 27 MAGs from the maize phylloplane (".B"), respectively.(b) Box plot showing the quality information (e.g., completeness and contamination) of all MAGs.(c) Heat map exhibiting the abundance of functional genes (based on KEGG Orthology) related to nitrogen cycling processes (e.g., nitrogen fixation, nitrification, and denitrification) in 58 MAGs.Detailed information on the 58 MAGs, including completeness, contamination, genome sizes, and taxonomic classification, can be found in the supplementary materials shared on Figshare (https://figshare.com/s/24f2a51568915f1138d6).(d) Nitrogen cycling pathways recovered based on the functional genes in 58 MAGs from the plastic and maize phylloplanes.The numbers in the brackets after the gene name represent the number of genomes identified and the total gene abundance, respectively.The number in the inner ring represents the average gene abundance (reads per million, RPM) based on the contig-level analysis.

FIG 6
FIG 6The conceptual models of crop diazotrophic community assembly.Crop-associated diazotrophic communities are largely shaped by environmental factors including climate (e.g., temperature and precipitation), edaphic factors (e.g., soil type, pH, and nutrients), and soil fertilization practices.These environmental factors define not only plant growth and development but also the microbial species pool and soil-plant interactions, which, in turn, shape crop diazotrophic community assembly.